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1. INTRODUCTION 



Concurrent developments in tropical meteorology, in particular the 
numerical simulation of tropical cyclones, fine-mesh modeling of extra- 
tropical systems and computer technology have made it feasible to con- 
sider a fully dynamic treatment of tropical cyclone motion and develop- 
ment. The purpose of this paper is to review these developments and 
discuss the feasibility of developing such a model for the western Paci- 
fic Ocean region. Development of an operational model for predicting 
landfall and dissipation of hurricanes over the Atlantic coast regions 
of the United States is proceeding at the Regional Modeling Branch of 
the National Meteorological Center. 

A smalT grid size is necessary to resolve the intense wind field and 
precipitation regimes near the center of typhoons. However, there seems 
to be evidence that the finest scales of motion may not have to be re- 
solved to get an acceptable forecast of tropical cyclone movement. It 
seems unlikely that the dynamic model will be able to improve on the 
short-term (less than 24 h) movement forecast techniques now in existence. 
This is primarily due to the inherent limitations in determining 
characteristics of each storm with present reconnaissance techniques. 

The dynamic model may be expected to be most useful in the cases in 
which the tropical cyclone is affected by extratropical flow. An example 
is the period of recurvature of the storm into the mid-latitude wester- 
lies, which presents one of the most severe problems for the typhoon 
forecaster. Consequently a dynamic model is likely to be advantageous 
in the 48-h to 72-h range. A possible exception to this may be the 
application of a dynamical model for sea and swell forecasts in the open 
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ocean, and for storm surge forecasts in coastal regions. Because the 
intense wind regions are highly concentrated, the inclusion of the 
atmospheric forcing to these models may require the space and time 
scales achievable only with a dynamical model. 

An operational numerical model capable of simulating the important 
features of tropical cyclones must have a much finer resolution than 
that required for extratropical cyclone prediction. As it is not 
presently feasible, or necessary, to maintain the same grid size far 
from the tropical cyclone, the general use of a grid(s) of finer reso- 
lution nested within the normal scale numerical weather prediction 
models is suggested. In the following sections, the status of nested- 
grid models and the problems involved in fine-mesh modeling will be 
reviewed. After a general discussion of the status of hurricane models, 
some consideration will be given to the feasibility of an operational 
model for the western Pacific Ocean. Such a developmental effort will 
require satisfaction of several interrelated tasks. These tasks will 
be described and alternative paths for attaining the objective will be 
suggested . 

2. STATUS OF NESTED GRID MODELS 

The need for finer scale models for tropical cyclone prediction has 
been described above. An awareness of the desirability of improved 
resolution in extratropical regions has also led to the development of 
regional models. In some cases, national needs may dictate a finer 
resolution than the usual 380 km resolution of a hemispheric model. 
Examples may be found in the efforts of several European countries to 
develop regional models suitable for their forecasting requirements. 
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More recently this has led to formation of a multi-national group to 
develop a fine-scale, regional model for the European area. 

There are both computational and physical reasons for developing 
these models. One may generally state that a mesh of 8-10 points is 
required to resolve adequately a meteorological feature, such as a 
travelling short wave. With fewer points, the phase speed will be 
decreased significantly and the amplitude of the feature may not be 
well-represented. This is particularly important in forecasting 
rapidly moving meteorological waves and the associated weather. The 
physical processes associated with small-synoptic scale and mesoscale 
atmospheric motions may also dictate a finer resolution. Although 
large-scale forcing associated with development of an extratropical 
cyclone leads to frontogenes is , the scale of motion of the front and 
its associated weather cannot be fully resolved on a synoptic-scale 
grid. Likewise, squall line development is triggered by larger scale 
features, but the processes involved in the generation and movement of 
the squall line are not represented by the normal governing equations 
for large-scale flow. In this regard the short-synoptic scale and 
mesoscale features present similar difficulties for numerical predic- 
tion as does the tropical cyclone. As these smaller scale forecast 
problems are common, the efforts of many groups have been brought to 
bear on the problems and have contributed additional expertise for solv- 
ing the tropical cyclone prediction problem. 

The primary feature which distinguishes nested grid models is the 
nature of the interaction between the grid scales. That is, how the 
information passes between two grids of different scales. A common 
technique, called one-way interaction, by Phillips and Shukla (1974), 
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is to first run the coarse-grid model saving the information on the bound- 
aries of the fine mesh. Then the fine-grid model is run using the coarse-mesh 
forecast as boundary conditions. It should be recalled that decreasing 
the grid size requires an associated decrease in time step. Thus, a space 
and time interpolation of the coarse-mesh values is required to provide 
the proper boundary conditions at some of the points along the fine mesh 
perimeter. It is clear that this process of reducing the grid size can be 
repeated until the desired resolution is achieved. Since the model is run 
independently for each grid size, and the same model code is normally used 
for each scale, this technique is easy to program and requires only the 
additional storage necessary to save the fine-mesh boundary values. 

There are two key limitations to this method. The inherent assumption 
in the one-way interaction models is that the larger scale motion determines 
the small scale motion, without significantly being affected by the processes 
occurring within the fine-mesh grid. Thus, the essential purpose of the 
finer resolution is to provide a better representation within a region from 
a purely computational viewpoint. An example is to decrease the phase error 
of small-synoptic scale waves. Another example of the technique is its use 
to provide synoptic-scale forcing for very fine resolution air-quality 
models, because the air-quality model is passive with respect to the synop- 
tic-scale model. Thus the effect of different synoptic-scale regimes can 
be imposed on the air-quality models, without the need for a feedback from 
the area with finer resolution. There is strong evidence that tropical 
cyclones, short-synopt ic and mesoscale features are not passive relative 
to the synoptic scales. These features trigger redistribution of energy 
and momentum that feeds back to the larger scale flow. Consequently it may 
be possible to forecast the initiation of these events with a one-way 
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interaction model, but forecasting the active phase of interaction with 
the larger scale is prohibited. In one-way interaction models the 
boundaries of the fine mesh must necessarily be far removed from active 
development regions. 

There is also a computational limitation which is common to all models 
of varying grid sizes, but it is particularly obvious for one-way inter- 
action models. The accuracy of finite-difference equations is dependent 
on the grid size and time interval. With a discrete change in grid size 
between two points, the solutions of the finite-difference equations will 
not be consistent. A wave that is well resolved on the fine-mesh grid, 
but not on the coarse grid, will have a different phase and amplitude as 
it progresses across a region. Consequently the boundary value provided 
by the coarse grid integration in a one-way interaction model must be 
incorrect. This is true for all inert io-gravity modes as well as for the 
meteorological modes that are being predicted. The result of the incon- 
sistent boundary conditions is a generation of spurious waves at all out- 
flow points of the fine-mesh grid. Since the proper mass, momentum and 
energy fluxes are not permitted at the boundary, the reflected waves will 
propagate into the interior and eventually destroy the solution. It might 
be noted that this would occur even if it were possible to provide the 
analytical solution to equations on the boundaries, as this would still be 
inconsistent with the solutions based on finite-difference equations. At 
any rate, some technique must be found to specify the outflow boundary 
conditions to prevent the reflection of the waves at the boundary. As will 
be described later, this is possible for some simple linear sets of equa- 
tions, but is generally not possible for all the wave modes permitted in 
the normal meteorological equations. 
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The two-way interaction models require a simultaneous integration on 
all grids. As the coarse grid integration provides the time-dependent 
values for the fine grid, and conversely the predicted fine-mesh values 
are a continually updated subset of coarse mesh, the solutions on the two 
grids are not permitted to develop separately. Because the coarse grid 
solution is influenced by the fine grid solution, the inevitable differ- 
ences at boundaries are controlled, without the heavy smoothing required 
in one-way interaction models. The need for simultaneous integration does 
require a more complex program than the one-way interaction technique, as 
the time step is different for each grid. With multiple reductions in 
grid size, and provisions for movement of the grids, a sufficiently fine 
resolution can be achieved. If the region of large gradients is relatively 
small, the horizontal array size in each grid does not have to be large. 
Because only the arrays for a single grid must be in the central core at 
a time, the code is efficient for a computer with a rapid transfer to 
external mass storage. This, in principle, is equivalent to solving the 
equations over only segments of the region. 

In summary, the one-way interaction technique can be used in situations 
in which the interior motions are passive relative to the large-scale forc- 
ing, or when the boundaries of the fine mesh are well-removed from the 
regions of redistribution of energy and momentum between the scales. A 
two-way interaction model is required where interaction between the differ- 
ent scales of motion occurs. If this region is localized, as in the 
tropical cyclone case, a movable fine grid with simultaneous integration 
should be used to properly resolve the scale interaction. Because an eight- 
fold increase in computer time is required to increase the horizontal reso- 
lution by a factor of two, it is essential that the fine mesh region be 
kept as small as possible. 
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One alternative to the discrete variation of grid size discussed above 
is the use of a continuously variable grid interval. Thus no abrupt 
change in phase and amplitude occurs across boundaries as in the discrete 
case. A variable resolution grid has been used in many (stationary) hurri- 
cane models to resolve adequately the inner region without having to compute 
at superfluous points in the outer region (e.g., Kurihara and Tuleya, 

1974). By specifying the functional dependence of grid size with radius, 
the equations may be transformed to an equivalent regular grid. The 
accuracy of the transformed equations varies with grid size. Furthermore 
the time step is that required for the smallest grid interval. However, 
the primary shortcoming of this type of model for operational use is the 
complexity for handling moving storms. As the center of the storm moves, 
the grid must also be moved to keep the large gradient region within the 
finest grid resolution. In practice, the grid may have to be moved nearly 
every time step, thus requiring interpolation in space and time of new 
boundary values around the perimeter of the region. Whether a realistic 
interaction between the hurricane circulation and the large-scale flow 
would occur through this interpolation has not been answered. 

The characteristics of limited-area on regional models can generally 
be grouped into four topics: (1) boundary conditions used for one-way or 
two-way interaction; (2) numerical considerations, including lattice arrange- 
ment of the variables and the time difference scheme; (3) initialization 
procedure, and (4) the representation of physical processes. A discussion 
of these topics with respect to a number of operational and research- 
oriented regional models is presented in the Appendix. These considera- 
tions will also form the basis of the discussion of three-dimensional 
tropical cyclone models in the next section. 
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3. STATUS OF NUMERICAL TROPICAL CYCLONE MODELS 



Two-dimensional models . There exists an extensive literature regard- 
ing the simulation of the symmetrical circulation in a hurricane. The 
motivation for these models is to simulate the important mechanisms that 
are involved in producing the intense winds near the center. The coupling 
between boundary layer convergence and associated latent heat release in 
deep convective towers near the center was shown to be an important aspect 
of hurricane development. A number of numerical problems, including the 
formulation of stable integration schemes, were examined in the process of 
developing these models. Although these two-dimensional models are not 
appropriate for real-time prediction, they continue to have great utility 
in testing new representations of physical processes, numerical methods, 
or the initialization of three-dimensional models (Anthes , 1974). It will 
be useful to review briefly the results of the two-dimensional models 
before proceeding to the three-dimensional models. 

The various simulations of hurricane development may be generally 
grouped into three phases — (a) organizational period, (b) rapid development, 
and (c) a quasi-steady state. The organizational period is necessary since 
the symmetrical models are started from a weak cyclonic system to avoid 
prejudicing the final state. Typically the maximum wind of 10 m/sec is 
located at several hundred km from the center. Because the flow is cyclonic, 
the boundary layer flow is directed toward the center. This leads to a 
boundary layer pumping as the air must ascend before reaching the center. 

The horizontal convergence of moisture flux plus the air-sea exchange 
would eventually lead to saturation of the boundary layer air, in the 
absence of convective fluxes through the top of the boundary layer. An 
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upward flux of moist air occurs within the clouds and drier air subsides 



in the environment. One of the important results of the symmetrical models 
is that convective release of latent heat occurring within the radius of 
maximum wind (RMW) tends to displace the RMW inward. Thus in the organi- 
zational stage the convective heat release increases the pressure gradient 
within the RMW by preferential heating of the columns near the center. As 
the pressure gradient is increased, the maximum wind increases. With a 
stronger cyclonic flow, the inward boundary layer flux strengthens. This 
provides additional fuel for convective heat release near the center, and 
the process continues. It is then clear that the rate of kinetic energy 
generation by the air flowing toward lower pressure must lead to develop- 
ment of a more intense system, unless the frictional dissipation is quite 
large. 

A rapid development phase occurs in all the symmetric models as the 
RMW approaches the center. Not only is the boundary layer inflow increas- 
ing as the cyclonic forcing increases, but the upward flux is occurring 
over a progressively smaller domain within the RMW. Very large heating of 
this narrow vertical column occurs, which rapidly deepens the central 
pressure, and increases the pressure gradient and wind speeds. The pro- 
cess of feedback between the convective heat release and the large-scale 
flow, which forces the boundary layer flow, is called Conditional Insta- 
bility of the Second Kind (CISK) . Thus, CISK guarantees development in 
these models (with realistic friction); however, the direct extension to 
nature may not be assumed since less than half of all potential disturb- 
ances develop into severe tropical cyclones. 

After the heating near the center brings the column to a moist 
adiabatic state, further deepening of the central pressure by convective 
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heat release is halted. Non-convective heating continues due to the 
larger scale ascent, but this simply maintains moist adiabatic ascent and 
does not warm the column. Thus the minimum surface pressure is attained 
and the heating moves outward into regions that remain convectively un- 
stable. If the model was run for sufficiently long periods the entire 
central region would become moist adiabatic, but the minimum pressure of 
the storm would remain in quasi-steady state. 

The importance of various initial conditions has been tested with 
symmetrical models (Rosenthal, 1971). A more moist initial state will 
shorten the organizational phase, but does not change the ultimate inten- 
sity of the storm. A more intense initial vortex will also enhance the 
development process in these models. The importance of latent and sensible 
heat has also been studied. Addition of heat and moisture from the ocean 
is crucial in raising the equivalent potential temperature of the inflowing 
air, which enhances the latent heat release that drives the storm. This 
is demonstrated in nature as hurricanes crossing over land lose this 
supply of moisture and eventually decrease in intensity when the primary 
moisture supply is only due to advection from outer regions that remain 
over water. As suggested above, the frictionally- induced inflow plays an 
important role in tropical cyclone development. Slight increases in the 
drag coefficient enhance the boundary layer pumping, and thus the convec- 
tive heat release, and the development time is shortened. However, the 
maximum intensity of the storm in the quasi-steady state is reduced by a 
larger drag coefficient, because the frictional dissipation is increased. 

In many of these respects, the symmetrical models provided guidance for 
the development of three-dimensional models that will be necessary for 
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intensity and movement forecasts of actual storms. It should be empha- 
sized that the two-dimensional models very realistically simulate the 
symmetr ical aspects of intense hurricanes. The complex wind and tempera- 
ture variations near the center are well-modeled, and thus one has some 
confidence that the mechanisms contained within these models are capable 
of reproducing these features in asymmetric models. However, the neces- 
sary resolution appears to be 10 km, if the intense wind region is to be 
resolved. This is crucial for three-dimensional models and leads to 
development of variable resolution, or nested-grid, models. 

Three-dimensional models . The first successful three-dimensional 
simulation of hurricane movement and development from real data was 
achieved by Miller (1969). This model was in many respects similar to that 
of Krishnamurti (1973), which has also been used to study non-developing 
tropical disturbances. Furthermore, a nested version of Krishnamurti ' s 
model has been developed by Mathur (1974), and used with real data to 
simulate hurricane development. The ability of this type of model to fore- 
cast developing and non-developing tropical disturbances represents an 
important achievement. Therefore it seems likely that this model contains 
enough of the essential physics to begin real-time forecasts, if the data 
and numerical problems can be overcome. The above models were initialized 
from hand-drawn charts and required extensive computer resources. Miller 
is continuing to develop a real-time forecast model (independent of the 
work at the Regional Modeling Group under J. Hovermale) using the NMC 
wind analyses from 850 to 100 mb as input to the model. 

The first fine-mesh (30 km) hurricane simulation capable of resolving 
the internal structure was described by Anthes , Rosenthal and Trout (1971). 
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Much of the physics of this model was also similar to the above models 
that used real data. Subsequent model development has continued under 
Rosenthal and R. Jones with the ultimate task of "simulating the track, 
intensity and inner core structure of tropical cyclones." In addition, 

"the model must also have the capability to provide rigorous simulation 
of various cloud seeding tactics associated with modification strategies 
proposed by Project Stormfury for hurricane moderation" (Rosenthal, Hur- 
ricane Research Progress Report #73, Nov. 1974). Finally, the fourth 
effort within NOAA to simulate hurricane development has recently been 
described by Kurihara and Tuleya (1974). This model incorporates a variable 
mesh grid ranging between 20 and 100 km resolution and includes an improved 
representation of physical processes. An important aspect of the improved 
physical aspects is that the model seems capable of simulating both the 
symmetrical component and the spiral bands. Because the vortex does not 
interact with a mean flow, this model is not capable of forecasting move- 
ment or development of actual storms. 

In this section the general characteristics of selected three-dimen- 
sional tropical cyclone models (see Table 1) will be described in the same 
manner as for the regional models discussed in the Appendix. As this 
material is primarily based on published articles, the information is 
somewhat dated. Nevertheless the general characteristics of the estab- 
lished models listed in Table 1 are probably not changing. The Regional 
Modeling Group (NOAA) model is currently being developed for real-time 
forecasts and is more subject to change. 

Boundary conditions . Results from the symmetrical models (e.g. 
Rosenthal, 1971) indicated that the boundaries should be well removed 
from the center of the hurricane to prevent contamination. Some recent 
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TABLE 1. Characteristics of selected three-dimensional tropical cyclone models. 



Principal 
Investigator (s) 
Group 


Nested? 

j. 


Boundary 

Conditions 


Number 
of grids 


Grid 
Size(s ) 
(km) 


#Time Dif- 
ferencing 


Vc >V 

Space Dif- 
ferencing 


Vertical 

Coordinate 


Layers 


Miller 

NHC-NOAA 


No 


C 


1 


75 


LF 


A 


P 


7 


Anthes & Rosenthal 
NHRL-NOAA 


No 


0 


i 


30 


EB 


A 


G 


3 


Jones & Rosenthal 
NHRL-NOAA 


Yes 


2 


3 


10/30/90 


EB 


B 


a 


3 


Harrison 

NPS 


Yes 


2 


3 


25/50/100 


LF 


A 


P 


3 


Mathur 

FSU 


Yes 


2 


2 


37/74 


QL 


A 


P 


3 


Madala 

FSU 


Yes 


2 


3 


20/60/180 


SI 


C 


z 


3 


Kurihara 

GFDL-NOAA 


No 


C 


1 


20/100 


EB 


A 


a 


11 


Hovermale 

RMG-NOAA 


Yes 


1 


2 


60/380 


LF 


A 


a 


10 


* 1 = one-way; 2 = 


two-way; 


CY = 


cyclic 


; 0 = open; 


C = closed 







# LF = leap-frog; LW = Lax-Wendroff ; SI = semi- implicit ; EB = Euler-backward 



**See text and Fig. 1 
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analytic results by Koss (Hurricane Research Progress Report #73, Nov. 
1974) show that the solutions are strongly dependent on lateral boundary 
conditions when a small domain is used. Koss suggests that a domain in 
excess of 3000 km may be necessary, based on his model solutions. How- 
ever, calculations with a fine mesh over such a large domain would 
require enormous computer resources. In order to proceed with actual 
forecasts, a much smaller domain has been used. Thus the lateral bound- 
ary conditions assume considerable importance, not only in representing 
the proper solution within the domain, but also for suppressing spurious 
computational modes. 

Miller's (1969) attempt to forecast hurricane development with con- 
stant geopotential values on the lateral boundaries resulted in unrealis- 
tic gradients of vorticity near the boundary. More recently Miller et al 
(1972) hold the u and v components of the wind constant on the two outer 
rows during the forecast cycle. After each step these two rows are ad- 
justed so the area means of divergence and vorticity for the entire grid 
are equal to those within the dynamically active regions. To constrain 
further the amplification of spurious waves near the boundaries, a large 
value of lateral viscosity is applied at the third and fourth rows in 
from the boundary. Open boundary conditions were essential for the small 
domain size (435 km radius) used by Anthes et al (1971) for simulation of 
hurricane development. Boundary conditions similar to those employed in 
symmetric models were used. Constant pressure and temperature values 
were specified around the boundary, and the horizontal velocity components 
were extrapolated to the boundary such that the horizontal divergence and 
relative vorticity tended to vanish at the boundary. Kurihara and Tuleya 
(1974) were able to use a 4000-km square domain because of the variable 



resolution grid. Consequently a closed boundary condition prohibiting 
any normal flux across the boundary was realistic. The normal component 
of the pressure gradient force was obtained from the geostrophic balance 
relation with the wind tangent to the wall. 

The only nested tropical cyclone model with one-way interaction listed 
in Table 1 is that of the RMG-NOAA. In the early stages of development, 
it was anticipated that this model would use boundary conditions similar 
to those of the limited fine-mesh (LFM) model at NOAA. That is, values 
along the boundary would be stored from the earlier run of the operational 
forecast model. If the storm was within the LFM domain the coarse grid 
size would be about 190 km at 60 N; otherwise the baroclinic primitive- 
equation model with resolution of 381 km at 60 N would be used. The 
stored values would be used to update the boundary values of the hurricane 
model, with perhaps some smoothing near the boundaries to damp the spurious 
modes. More recently this group has been testing two-way interaction 
models of the type to be described below. 

One-way interaction boundary conditions must also be specified for 
the coarse-mesh grid (CMG) of the nested models. Mathur (1974) held the 
initially analyzed variables constant throughout the forecast period, 
and applied a nine-point smoother on the grid points within three rows 
from the boundary. The particular situation that Mathur considered in- 
cluded little interaction between the developing hurricane and the 
surroundings. In cases with stronger interaction, the constant boundary 
conditions would have been inappropriate. As Harrison (1973) used 
analytical initial conditions, the CMG boundary conditions were chosen 
appropriate for a channel flow with cyclic conditions at the east-west 
boundaries. Allowing no flow across the northern and southern walls 
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tends to emphasize the westerly movement of the disturbance. Madala 
(1973) also used a channel flow for the coarse grid, but he used up- 
stream differencing to obtain the necessary boundary values. 

The key aspect of the nested grid integrations is the two-way in- 
teraction across each discrete change in grid size. Harrison (1973) 
adapted the one-dimensional method of Harrison and Elsberry (1972) for 
specifying the variables on a two-dimensional FMG from the computed 
CMG tendencies along the interface. This method was successful in simu- 
lating the intensification of a weak tropical circulation, and in parti- 
cular, allowed the transfer of mass and energy through the interfaces 
as the system developed. However, the method has not been tested with 
the more intense stages. Jones and Rosenthal (personal communication) 
have used the Harrison method in successful simulations of hurricane 
development. Mathur (1974) also used the tendencies at CMG points on 
the interface common to both the FMG and CMG in calculating the pressure 
gradients and nine-point, quasi-Lagrangian advection terms. Rather than 
interpolating the CMG tendencies in space to the non-coincident FMG grid 
points along the interface, as in Harrison (1973), the variables were 
assumed to have continuity in the second-order derivatives calculated 
at the adjacent CMG points. The latter step tends to smooth the values 
along the interface between grids. However, since all variables are in- 
terpolated in this manner, the variables are not necessarily in dynamic 
balance at the intermediate points. 

In the semi-implicit model of Madala (1973), a governing Helmholtz 
equation for the geopotential field must be solved for each grid (see 
Appendix). Madala simultaneously solved the Helmholtz equations for the 
three grids. Furthermore, new boundary conditions for the inner grids 
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were calculated (from the coarse-grid solution using quasi-Lagrangian 
advection) after each iteration of the relaxation process. Presumably 
the solution would have converged more rapidly if the coarse-grid solu- 
tion would have been calculated, and boundary conditions set only once 
for each successively smaller grid. After the finer mesh solution is 
obtained, the coincident points in the interior of the CM G would be re- 
placed by an average of the surrounding points in the FMG. This would 
produce an interaction directed from the FMG to the CMG comparable to 
the explicit solution method of Harrison (1973). Interaction between 
the different grids occurs in every iteration of the Helmholtz equation 
solution in Madala's case, rather than during each time step as in 
Harrison's method. This highly interactive solution method may account 
for Madala having to reduce the time step by a factor of more than 10 
during the course of the integration. It is not clear whether the 
method suggested above with interaction only during a time step would 
give a comparable solution. If so, this procedure would be more 
economical, as a longer time step could be utilized in each grid 
interval . 

In summary, the boundary conditions for a model capable of resolving 
the inner features of a tropical cyclone must provide for two-way inter- 
action in a nested-type model. This allows mass, momentum and energy 
transfer and reduces computer requirements compared with a uniform 
fine-mesh region over the entire domain. 

Numerical considerations . The time differencing schemes for the 
models listed in Table 1 are as varied as in the fine-mesh, mid-latitude 
models (see Appendix). Leapfrog differences are used for simplicity 
in the Miller et al (1972) and Harrison (1973) models. The NHRL and 



17 



GFDL models use Euler-backward differencing to damp the high-frequency 
modes. Mathur (1974) has used the quas i-Lagrangian scheme for evaluat- 
ing the substantial derivatives. This scheme has advantages in minimiz- 
ing phase errors over the normal centered-space and time technique. 

The importance of this advantage is not obvious since the reason for 
using a nested grid is to provide adequate resolution for the short 
waves. On the other hand, the grid size may not have to be reduced as 
much to obtain the same phase error properties if the quas i-Lagrangian 
scheme was incorporated. Furthermore, use of the quas i-Lagrangian (or 
higher order space differencing) would be useful in the coarse-mesh 
grid to provide optimum handling of the phase of the large scale waves. 

An efficient, semi- implicit time differencing scheme has been 
developed by Madala (1973). The pressure gradient and divergence terms 
are treated implicitly, because these terms control the high-speed 
gravity waves. If the time necessary to solve the Helmholtz equation 
could be further reduced by newly developed techniques, the advantages 
of the semi- implicit technique could be exploited. In some cases these 
new techniques have constraints on the number of grid intervals, 
although this should have little impact on the development of a new 
model. 

A majority of the nested tropical cyclone models listed in Table 1 
carry all variables at each gridpoint (Scheme A in Fig. 1). Although con- 
venient, this scheme is inefficient and results in noise generation 
as the solutions on different lattice points tend to separate (Gerrity 
and McPherson, 1970). Madala (1973) used Scheme C, which Arakawa et al 
(1974) suggest has favorable properties both for geostrophic adjustment 
and for simulating the slowly changing quasi-geos trophic modes. Madala 
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notes that this scheme is also very desirable in formulating the Helm- 
holtz equation that arises with semi-implicit time differencing. 

Elvius and Sundstrom (1973) showed that an even more compact form 
which required fewer Helmholtz equation solutions was possible. Thus, 
the interrelation between optimum space and time differencing should 
be considered in planning a fine-mesh tropical cyclone model. 

For open-ocean forecasting of tropical cyclones, it would seem 
appropriate, as well as convenient, to use a pressure-coordinate model. 
However, forecasting landfall and the subsequent filling stages of 
tropical cyclones may require the sigma coordinate system. This is 
probably most important if a prime consideration is for precipitation 
forecasts as the storm moves over mountainous terrain. Furthermore, 
if the hurricane model boundary conditions are to be derived from a 
hemispheric model in sigma (a) coordinates, it seems essential that 
the hurricane model be cast in a coordinates. Thus the use of vertical 
height (Z) coordinates by Madala can be viewed as a disadvantage for 
initialization and for operational, nested-grid forecasts over terrain. 
A a coordinate version of Madala' s model would appear to be most 
favorable for tropical cyclone forecasting. 

Ini t ia 1 iz at ion . The discussion of initialization for fine-mesh, 
mid-latitude models (see Appendix) is also applicable to the tropical 
cyclone problem. Certainly the emphasis must be upon an accurate 
specification of the wind field. However, the scarcity of data, parti- 
cularly on the scale of intense wind regions, will require an optimum 
use of all available data. It remains to be seen if useful information 
can be obtained from satellite data. The tropical storm is typically 
covered by a thick overcast so that sensors in the carbon dioxide bands 
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cannot be used to determine the vertical temperature profiles. Another 
difficulty with polar-orbiting satellites is that the data are asynoptic 
and methods must be developed to blend this information to produce an 
initial field. If sensors (or processing techniques) could be perfected 
that would permit inference of the profiles on the desired scale, say 
20-30 km, the temperature errors on these scales may well be of the same 
order as the gradients. Even in the most optimum case, one must cor- 
rectly specify the wind field that is consistent with the deduced mass 
field. For once the integration is begun, the mass field will tend to 
adjust to the specified wind field. 

Given the need for accurate wind information within and surrounding 
the tropical cyclone, maintaining the existing reconnaissance is an 
absolute minimum. The technology for rapid transmission of this recon- 
naissance data to the operational centers must be pursued. However, 
aircraft can only provide information within 150 km of the center at 
one or two levels. Perhaps extensive exploitation of satellite-derived 
cloud motions will provide low-level and upper-level wind information 
to mesh with the aircraft-sensed region. This near-environment wind 
regime will be essential for proper short-term steering of the model 
storm. Aircraft reports play an important role in tropical wind field 
analysis, but these tend to be concentrated near 250 mb. Between this 
layer and the surface, where ship reports tend to improve the coverage, 
there are few wind reports. 

The successful hurricane forecasts with real data (Miller et al, 
1972; Mathur, 1974) have used hand-analyzed wind fields. This has the 
advantage that the analysts subjectively interpolated into data gap 
regions to provide a realistic analysis. An important point is that 
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Mathur (1974) analyzed on the fine mesh scale throughout the region and 
discarded the unnecessary points within the coarse grid. Miller (1969) 
also had initial data on the fine scale. The alternative, interpolation 
of large-scale data to smaller scales, has not been tested with real 
data cases. Clearly some additional information must be provided for 
an accurate representation near the center, but the amount and distribu- 
tion are not obvious. It should also be noted that Mathur' s and Miller's 
data were analyzed at each vertical level, and presumably were carefully 
checked for vertical consistency. As indicated above, the available 
data sources tend to be concentrated at upper levels and the surface. 

Some attention must be given to producing objective analyses that are 
vertically consistent between these levels. 

As Anthes (1974) points out, because the local forcing is so strong 
in a tropical cyclone, there is hope that the initial data need not be 
of high resolution over the entire horizontal and vertical extent of the 
domain. Anthes used a symmetric model of a nearly steady-state hurri- 
cane to determine the variables that are most important in providing an 
adequate initial state. The low- and middle-level observations were most 
useful in reducing the initial errors in Anthes' experiment. As would 
be expected, observations in the region of highest winds were crucial. 

If only observations from the surroundings were utilized, the capability 
of the model to reproduce the control state was rather poor. Note that 
these variables may not be the same as in the situation in which the 
development of the hurricane is to be simulated, rather than the dupli- 
cation of an existing intense hurricane. That is, a more (or less) 
stringent specification of initial conditions may be applicable depend- 
ing on the purpose of the forecast. 
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Miller et al (1972) normally start their model with wind fields. 

The height fields can be determined from a version of the divergence 
equation, and temperatures can be computed from the thickness fields 
using the hydrostatic equation. Vertical motions are diagnosed from 
an omega equation and the divergent component of the wind field is 
derived. The temperatures and winds are then forecast. After a short 
period, the rotational component of the forecast wind fields is retained 
as initial data, and the initialization cycle is repeated to obtain a 
new mass, temperature and divergent wind component fields. This diag- 
nostic-prognostic cycle is repeated four or five times until the 
heights and vertical motions are relatively stable. Mathur (1974) 
followed a similar procedure, but did not repeat the initialization 
cycle. Mathur used a balance equation and a quas i-geos trophic o> equa- 
tion that are subsets of the complete divergence and cu equations of 
Miller et al. Most of the models listed in Table 1 were started from 
a weak vortex that was in gradient balance, and the initialization 
phase was of little interest. 

Because the amount of data available within the tropical cyclone 
region will normally be insufficient, some method is necessary to in- 
troduce, or bogus, data on the proper scales. For example, the actual 
values of maximum wind and the eyewall radius available from satellite 
or reconnaissance data are insufficient, but one may be able to con- 
struct a model that is consistent with these storm parameters. At the 
RMG, a two-dimensional version of the complete model will be run to 
produce initial conditions that are consistent with the environmental 
sounding. Such a symmetrical model would require little computer time 
to achieve steady conditions, and the symmetrical component of the storm 
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can then be used as bogus data during the initialization phase. In a 
four-dimensional assimilation scheme, the complete model equations are 
integrated from a prior (say 12-h old) time. During this initialization 
integration, a correction will be added at each time step to "nudge" the 
predicted gridpoint value toward the "observations . " The weighting factor 
for the correction term may be a variable which "decreases with observa- 
tion error, increasing horizontal and vertical distance separation, and 
increasing time separation" (Anthes, 1974). The advantage of this system 
is that the observations never replace the forecast variables and are 
never added impulsively. In some experiments with a constant weighting 
factor, Anthes showed that only a few selected points of data were suffi- 
cient to dynamically initialize a symmetrical hurricane model. 

It may also be possible to bogus data derived from empirical models 
of the symmetric component in mature tropical cyclones. For example, 

Riehl (1963), LaSeur (1966), Shea and Gray (1973), and Gray and Shea 
(1973) summarize observations of the two-dimensional structure of hurri- 
canes. Given only a few characteristics of the actual storm, one could 
specify an inner structure that is consistent with these observational 
models. Elsberry, Pearson and Corgnati (1974) have used such a quasi- 
empirical model to determine the boundary layer characteristics and the 
maximum intensity of the model storm for different sea-surface tempera- 
tures. This type of model could be used to specify the symmetric wind 
and thermal structure that would be consistent with actual observations. 
Selected data from the model could then be used to enhance the data 
coverage for statically or dynamically initialized models. 

As indicated in the discussion of initialization for regional models 
(see Appendix), the variational approach of Sasaki (1958) can make use 
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of several different types of data, and blend the information into a 
form that is dynamically consistent with the model equations. These 
are desirable properties for an operational tropical cyclone model, as 
such a method would reduce the computer time necessary for dynamic 
initialization techniques. As the final state is not completely in 
agreement with the dynamical equations (or the observations), an adjust- 
ment will occur during the early stages. The effects of this adjustment 
can be reduced by properly specifying the weighting factors on the 
dynamical terms. 

Representation of the physical processes . The primary process re- 
quired in a simulation of tropical cyclone development is the latent 
heat release, which is related to the convergence of moisture within 
the boundary layer. Some crude representations of this highly compli- 
cated interaction have been used to demonstrate the role of latent heat 
release. Nested grid models by Harrison (1973) and Madala (1973) simply 
specified the distribution of latent heat release in the vertical as a 
function of boundary layer convergence. These models were not intended 
to duplicate observed cases of hurricane development, but this represen- 
tation was sufficient to demonstrate that nested grid models could 
properly handle cases of energy redistribution from the fine mesh to the 
coarse mesh. Progressively more complex representations are available 
for latent heat release and planetary boundary layer effects. Miller 
et al (1972), Mathur (1974) and Anthes et al (1971) have used modified 
versions of the Kuo (1965) parameterization of convective latent heat 
release. Kuo related the amount of heating to the ratio of the moisture 
supply to the moisture required to saturate the atmospheric column and 
bring it to the temperature of the moist adiabat passing through the base 
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of the cloud. The supply of moisture to the clouds was related to the 
net convergence of moisture into the column by the large-scale flow. 

In a tropical cyclone this convergence is primarily due to boundary layer 
convergence plus evaporation. The second key element in any convective 
parameterization scheme is the specification of the vertical distribution 
of the heat release. Kuo related the heating to the temperature differ- 
ence between the environment and the moist adiabat, which leads to a 
maximum heating rate in the upper troposphere, as is observed in actual 
storms. This representation of convective heat release seems to contain 
the necessary ingredients for hurricane development forecasts in both 
two- and three-dimensional models. The ultimate intensity of the model 
storm is the dependent on the moist adiabat near the center, because the 
Kuo heating function is self-limiting for a moist adiabatic environment. 
Continued maintenance of the storm requires a non- convective heat release 
proportional to the large-scale ascent rate. It is well-known that 
unconditional heat release will result in grid scale overturning, so 
this process must be controlled by smoothing (Miller et al, 1972). 

The most sophisticated representation of physical processes is found 
in the recent model of Kurihara and Tuleya (1974). Several aspects of 
the model are similar to the Geophysical Fluid Dynamics Laboratory general 
circulation model, although radiative effects are neglected. Subgrid- 
scale vertical diffusion processes include an air-sea interaction within 
a Monin-Obukhov framework, free dry convection and mechanical turbulence 
above the surface layer. Free moist convection is included if the 
vertical change of temperature in a rising cloud element is greater than 
the lapse rate of the cloud environment. The fractional entrainment rate 
for the cloud elements is assumed to be inversely proportional to the 



25 



radius of the cloud (only one size is permitted). A key assumption 
appears to be a specification of the efficiency of the convective trans- 
port by an empirical coefficient. The parameterization scheme is 
similar to that of Ooyama (1971), in that an arbitrary factor determines 
the importance of each cloud size. A closed solution to the parameteri- 
zation problem has been proposed by Arakawa and Schubert (1974). This 
new scheme will be tested by Jones and Rosenthal (Hurricane Research 
Progress Report, No. 73), since the tests of Project Stormfury hypotheses 
will be more convincing if some of the arbitrariness of the convective 
heat release rates is removed. Schubert (personal communication) also 
will test the new parameterization scheme for prediction of various tropi- 
cal weather systems. 

The progress made in representing the cumulus heat release, and the 
previous applications to hurricane development simulation, indicate that 
methods for representing these physical processes are sufficient for an 
operational model. The degree of sophistication necessary, of course, 
depends on the purpose of the model. In the initial stages of model 
development, the Kuo scheme would seem adequate, as considerable exper- 
ience has been gained with this method of parameterization. A compari- 
son with more sophisticated parameterization schemes would be appropriate 
at a later stage. 

4. PROPOSAL FOR A TROPICAL CYCLONE MODEL DEVELOPMENT 

As noted in the above review, there is active research in the area 
of tropical cyclone simulation. We propose development of an operational, 
fine-mesh tropical cyclone model suitable for Navy needs in the western 
Pacific Ocean. The considerable progress toward solution of the numerical 
aspects of this problem puts this research in the category of research 
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applied to Navy needs. It is recommended that this effort make maximum 
use of concurrent efforts among government (and non- g over nmen t ) agencies. 
Whereas the problems associated with the development of such a model are 
similar, and the basic model will have many common aspects, the required 
resolution, initial data and forecast interval will be determined by the 
different tasks. Exchanges of information with the various government 
agencies active in tropical cyclone prediction would assist in the develop- 
ment of a Navy model for WESTPAC (see Table 2). It seems essential that 
the proposed development be closely coordinated between the Environmental 
Prediction Research Facility, Fleet Numerical Weather Central and the 
Naval Research Laboratory. Furthermore, the sustained interest, and 
operational data inputs of the Fleet Weather Central/Joint Typhoon Warning 
Center will be necessary throughout the developmental period. The timing 
of this development effort, and the eventual operational application, will 
depend on the resources available to each group. In many cases the efforts 
dedicated to the tropical cyclone problem will be transferable to other 
regional, or global models; for example, any improvement in numerical 
techniques or representation of the physical processes. Exploiting the 
commonality then requires an awareness of the similar problems, which can 
be achieved through good communication, rather than a redirection of effort. 

The tasks to be undertaken to develop an operational model are listed 
in Table 3. For convenience, the developmental progress is grouped into 
four phases, although there inevitably will be some overlapping. Time 
sequences have not yet been established, because the level of funding is 
uncertain (a fully committed effort, which would include adequate computer 
resources for developmental and quasi-operat ional phases, perhaps would 
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TABLE 2. Possible coordination with other government agencies. 



Group 

NHC-NOAA 


Exchange of data/model comparison 


RMG-NOAA 


Exchange of data, model development and initialization techniques 


NHRL-NOAA 


Model development comparison/parameterization techniques 


NRL 


Model development for semi-implicit technique 


FWC/JTWC 


Data collection/evaluation 

Possible coordinated EPRF, FNWC and NPS endeavors 


EPRF 


Development of semi-implicit version with efficient Helmholtz 
equation solution. 

Initialization with variational techniques* 

Evaluation of possible satellite data input for WPAC. 
Parameterization of latent heat release development /evaluation. 


FNWC 


Latent heat release parameterization package from global 
model- evaluation. 

Planetary boundary layer package from global model. 

Sigma coordinate model 

Objective analysis package 

Initialization with variational techniques 


NPS 


Consideration of numerical aspects including time and space- 
differencing schemes. 

Evaluation of use of model data for bogusing. 

Evaluation of data and model requirements. 

Dynamical initialization experiments. 



Coordination with university or non-U. S. groups might also be considered. 
Groups at Colorado State University and Pennsylvania State University are 
active in fine-mesh modeling of tropical disturbances, and the Japanese 
Meteorological Agency has developed an operational, nested grid model. 
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TABLE 3. Proposed plan for development of an operational tropical cyclone model 
for the western Pacific region. 

Feasibility Development Development Quasi-operational 



Task 


phase 


Phase I 


Phase II 


testing phase 


Initialization 


Hand-analyzed 
wind data from 
WPAC typhoons . 


Hand-analyzed 
fine mesh data 
variational 
analysis 
technique . 


Objective analysis 
data tests. 


Improvement of 
ob j ect ive analys is . 




Test inclusion 
of divergent 
component . 


Bogus with 
dynamic / 
empirical 
model data. 


Comparison of 
dynamic and varia- 
tional techniques . 


Evaluation of 
bogus ing techniques 






Interchange 
of data with 
other groups 


Satellite data 
input. 




Numerical 

aspects 


Averaged 
pressure- 
gradient to 
reduce time. 


Staggered grid 
to reduce core 
requirements . 


Increase number 
of layers. 


Reduction in grid 
size to 20/60/180 
km. 




Triply nested 
model with 
grid sizes 
55/110/220 km. 


Semi-implicit 
time differenc- 
ing. 

a coordinate. 


Mesh with FNWC 
global model. 




Physical 

processes 


Simplified 

heating/BL 

package. 


Modified Kuo 

heating 

package. 


Improved heating 
package including 
radiation. 


Relate to oceano- 
graphic products. 






Improved BL 
representation. 






Goals 


Demonstrate 
feasibility 
for 36-h 
motion fore- 
cast of tropi- 
cal cyclones. 


Extend move- 
ment forecast 
to 48h 

for tropical 
cyclones . 


Extend movement 
forecasts to 72h 
for tropical 
cyclones . 


Evaluate intensity 
and movement fore- 
casts to 72h. 




Development of 
test model. 

Collect sample 
of test cases . 


Improve effi- 
ciency of 
test model. 

Evaluate data 
and model 
requirements 
through inter- 
change of 
test cases. 


Evaluate model 
forecasts for 
tropical 
disturbances . 


Near real-time 
forecasts in 
support of Project 
Stormfury and 
comparison with 
NHRL model. 
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require a three-year period to fully accomplish the goals listed at the 
bottom of the table) . 

The first phase may be thought of as a preliminary effort to demon- 
strate the feasibility of using a currently available model (Harrison, 
1973) to forecast typhoon movement* Data for one experiment have been 
prepared by FWC (GUAM) JTWC personnel. The meteorological situation is 
a late-season typhoon (IRMA), that did not recurve as expected. Such a 
potential recurvature situation is considerably more difficult than the 

case of Mathur (1974). As in Mathur's case, the initialization will be 
from the wind field. One purpose of this phase is to evaluate the ade- 
quacy of the data sources in the Western Pacific Ocean for specifying 
the environmental conditions surrounding a typhoon. As the case has been 
subjectively analyzed by an experienced forecaster, one might consider 
that this is an optimum analysis, at least compared to an objective 
analysis case. A second purpose will be to test the role of a fine- 
mesh analysis (or lack thereof) in forecasting movement. A procedure 
for establishing the initial center location within the grid will be 
required. A final goal is to establish a test model and procedure as 
a reference base for evaluting future improvements or modifications. 

The basic model of Harrison (1973) will be modified by time averaging 
of the pressure-gradient term (Brown and Campana, NMC technical note). 
This should permit a doubling of the time step, and thus allow more 
cases to be efficiently run. 

Based on Mathur's one case, and the work of Miller et al (1972) with 
a non-nest ed model, it is expected that the first phase will indicate 
the potential of a dynamic model forecast. One cannot satisfactorily 
demonstrate feasibility with a single case. Nor are we anticipating that 
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the Harrison model will be adequate as an operational model, since more 
efficient time and space differencing schemes are now available. Further- 
more, the model includes minimal physical effects. 

The primary efforts in developmental phase I should be directed to- 
ward initialization and numerical aspects of the model. Certainly addi- 
tional cases should be prepared by FWC/JTWC personnel. However, it is 
quite expensive to prepare even a single case, and a parallel effort 
should be made to establish a data exchange with other interested groups. 
As these groups have predominately used Atlantic cases, it might be 
beneficial for them to have access to additional cases from the Pacific. 
This exchange should also assist in comparisons of initialization tech- 
niques, and of the individual models. A considerable effort should be 
made to develop initialization techniques that are appropriate to the 
tropical cyclone case. The developing expertise in both variational and 
dynamical initialization methods should be utilized. It may occur that 
these more complicated techniques will not be required for satisfactory 
movement forecasts, although the intensity forecasts will likely require 
a superior initialization technique. The second emphasis in phase I 
should be the improvement of the model by exploiting recent developments 
in numerical techniques. Here we hope to make use of the experience 
gained by various active groups. For example, the NHRL group under Jones 
and Rosenthal are tasked with developing a model in support of the 
Project Stormfury experiment scheduled for the western Pacific area 
during 1977-78. As their primary goal will be to test typhoon modifica- 
tion hypotheses, it would be an advantage to them to have a comparable 
numerical model for 48- to 72-h forecasts of typhoon movement and inten- 
sity. This would greatly aid in their operational planning, and would 
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provide a reference to evaluate the seeded versus non-seeded experiments. 

A number of other possibilities for coordination with other groups are 
suggested in Table 2. By the end of the development phase I, enough test 
cases should be run to demonstrate the possibilities for an operational 
model-both in terms of data and an efficient model code. 

During development phase II the main objective may be described as 
a refinement of all aspects of the planned model. Evaluation of the 
initialization task should include a comparison of an objective analysis 

rk 

versus a subjective analysis. Within this time frame the potential use 
of satellite data over the western Pacific region should be pursued. 

The refinement of the model should be directed toward making it compatible 
with the operational FNWC system, so far as is possible considering the 
different goals of the tropical cyclone system. An increase in computer 
resources should be planned for this period as the number of levels is 
increased and the representation of physical processes is improved. 
Furthermore, the goal of extending the movement forecasts to 72 h will 
require additional computer time. 

Movement to a quas i-operational phase will primarily be dependent on 
new computer resources. The additional goal of forecasting intensity 
as well as movement will require a finer resolution. Essentially another 
order of magnitude in computer power will be required for the planned 
reduction in grid increment. Without the additional resources, one might 
still opt for movement forecasts, and perhaps some definition of the 



*This is not to say that a project designed to compare objective and 
subjective analyses could not be initiated during phase I as a number of 
test cases become available. 
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envelope of gale-force winds. It might be noted that a reduction of 
grid size would likely improve the movement forecasts as well* To obtain 
the full benefit of the dynamic model for fleet use, the output should also 
be tailored to the FNWC oceanographic products. Because of the United 
States Navy meteorological and oceanographic needs in the western 
Pacific area, it is hoped that the Navy community will give full support 
toward development of an operational, tropical cyclone forecast model. 
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APPENDIX 



PROBLEMS IN FINE-MESH MODELING 

Boundary conditions . The general problem of selecting boundary condi- 
tions for a regional model has been studied for many systems of meteoro- 
logical equations (e.g. see summary in Harrison and Elsberry, 1972; and 
recent work by Elvius and Sundstrom, 1973). For nested models, the selec- 
tion of a one-way or a two-way interaction model is dependent on the 
meteorological parameters to be forecast. If one develops a regional model 
independent of a hemispheric or global model, the boundary values must be 
generated from the internal flow or be held fixed. The early version of 
the limited fine mesh (LFM) model of the National Meteorological Center 
used constant boundary values (Howcroft, personal communication). Platzman 
(1954) showed that constant boundary values permit a computationally stable 
interior solution, but the noise generated at the outflow points propagates 
into the field and the solution is useless after a relatively short time. 
Heavy smoothing is normally required near the boundaries to keep the pre- 
dicted solution consistent with the constant boundary values. Krishnamurti 
(1973) uses constant boundary values after a 24 h adjustment period, during 
which a condition for controlling the horizontal divergence is applied at 
the first interior gridpoints. Thus the normal (meridional) component at 
the northern wall (J=M) may be written from the horizontal divergence 
expression 

v (I,M) = n[v(I,M-2) - ^|u(I+l,M-l) - u(I-l,M-l)J ] 

24-t 

where r\ “ ~ for t < 24 h and r\ = 0 for t > 24 h. A similar ex- 
pression is used for the meridional component along the southern wall. 
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TABLE 4. Characteristics of selected limited-area models. 
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The remaining variables along the walls are held fixed at the initially 
analyzed values. Since the Krishnamurti (1973) model was developed for 
short-term prediction in the tropics, the east-west boundary condition 
was made cyclic by inserting a buffer region across which the eastern and 
western boundary values are matched. Krishnamurti has used the foregoing 
boundary conditions to allow for propagation of inertio-gravity oscilla- 
tions at the boundaries during the 24 h adjustment phase. After 24 h 
the interior flow has adjusted to the presence of the walls. 

If the boundary values along the perimeter of the fine mesh are avail- 
able from a coarse-grid model, the interior solution will not be consistent 
with the boundary values. One may think of the problem as an over-speci- 
fication of the fine-mesh boundary values, which leads to a spurious 
computational solution. Charney (1962) has shown that for a barotropic 
equation set, the normal velocity component should be specified at all 
boundaries, and the potential vorticity at the inflow points. At outflow 
points, the vorticity may be extrapolated from the interior. For the non- 
linear barotropic model (shallow-water equations), both Elvius and Sunds- 
trom (1973) and Chen and Miyakoda (1974) propose specifying the normal 
velocity at all points and the tangential velocity at the inflow points. 

The remaining required variables at the outflow points are determined from 
the governing interior equations using one-sided differences, or by a 
method of characteristics (Chen and Miyakoda). These boundary conditions 
have been shown to be stable for the shallow water equations, but will be 
ill conditioned with a baroclinic model using the hydrostatic approximation 
(Oliger, personal communication). Thus in most cases of one-way interaction 
boundary conditions, an over-specification is used and various artifices 
are employed to control the resulting computational modes * In the British 
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model (Jones, 1974) the height changes (and other variables) on the fine 
mesh are saved from the earlier octagon model run. However, the normal 
velocity changes on the boundary are replaced by geostrophically calculated 
values. Strong diffusion is applied to the three rings within the boundary. 
Williamson and Browning (1973) describe boundary conditions for a limited 
area model (LAM), which may be located at any point within the National 
Center for Atmospheric Research (NCAR) general circulation model. If the 
wind computed by the global model (stored on tape) is directed into the LAM, 
the variables from the global model are used on the LAM boundaries. If the 
global wind is outward, an extrapolation scheme from within the LAM is used 
to determine the LAM boundary values. The extrapolation technique follows 
Shapiro and O f Brien (1970) in that a path, or space-time line, is calculated 
to determine the interior properties that will arrive at a boundary point 
during the time step. That is, the global winds are used to advect the LAM 
variables to the boundary grid point. Another technique to melt the 
interior solution to the external solution is to use so-called "spongy 
boundary conditions" (see Kesel and Winninghoff, 1972). For example, Kreitz- 
berg et al (1974) applied these conditions to the predicted tendencies of 
the property A, 



SA(B-n) SA(B-n) 

St “ 7(n) St 

where 7(n) = 0.3, 0.7 and 0.9 at n = 0, 1, and 2 points from the boundary 
point (B). These conditions are designed to control reflections at the 
boundary by permitting only a slow rate of change near the boundary. Never- 
theless, Kreitzberg et al (1974) found it necessary to apply a smoother- 
desmoother to the variables near the boundaries as frequently as every 5 
time steps. The "sponge" can also be used to provide a smooth transition 
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between two scales by the relation 



A j z ( B " n ) = 7 (n) A (B-n) + [l-7(n)J A (B-n) 

actj line coarse 

to all the variables. 

In summary, a number of techniques are available for using coarse- 
mesh values at fine-mesh boundary points in one-way interaction type 
models. In every case the spurious computational models generated near 
the boundaries must be controlled or the interior solution will be con- 
taminated. The smoothing, or filtering, must be carefully applied to 
avoid destroying the meteorological modes as well. If the coarse-mesh 
grid (CMG) is sufficiently large relative to the tropical cyclone, the 
CMG may be a subset of the hemispheric/global model, and one-way inter- 
action boundary conditions cf the tyre described above will be useful. 

In a two-way interaction model the time change at the outflow (and 
inflow) boundary points is calculated from the numerical solution at the 
corresponding CMG point. Consider the cyclic grid network depicted in 
Fig. 2 which changes spatial increment at the points j = 4 and j = 10. 

The ratio of the spatial increment in the CMG to that of the FMG is two, 
so the even numbered gridpoints for A < j < 10 are also coarse mesh grid- 
points. Thus a centered space difference across point 4 will involve 
points 3 and 6, although the solution at 6 is calculated on the FMG. The 
adjacent interior point j = 5 then uses the value at j = 4 calculated from 
the CMG solution as a boundary value in the centered space difference on 
the FMG. Two-way interactions specify CMG values as boundary conditions 
as in the one-way type, but the CMG solution is simultaneously calculated 
and interacts with the FMG through use of the interior values j = 6 and 
j = 8. Because the time step for the FMG is half that of the CMG, the 
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time sequence of the two solutions must be maintained with two steps in 
the FMG for every CMG time step. In some cases the same time step has 
been used in both grids, and although convenient for programming, this 
procedure is quite wasteful if more than two nested grids are used. Ex- 
tension of the above procedure to two dimensions may be examined with the 
aid of Fig. 3, which illustrates a multiple reduction of grid scales. 

Again the ratio of grid interval between successive grids is a factor of 
two (for staggered grids, the ratio is normally three). If each of the 
three grids shown contained the same number of gridpoints, the medium-mesh 
grid (MMG) would cover 1/4 of the CMG area, and the FMG only 1/16 of the 
CMG area. Thus the advantage of the nested grid is obvious if the impor- 
tant features to be forecast are contained within a small region, as in 
the case of a tropical cyclone. 

Numerical considerations . One of the common problems of every opera- 
tional forecast center is the need to generate the best forecast in the 
shortest time possible. Because each reduction of grid interval by a 
factor of two requires an eight- fold increase in computing time, the use 
of fine-mesh grid models has promoted extensive research into increasing 
the maximum time interval for a given grid size (the reader is referred to 
numerical weather prediction texts, e.g. Haltiner 1971, for discussion of 
the criteria for various time-differencing schemes). Several of the models 
listed in Table 4 use the standard, centered time difference, or leapfrog 
(LF) scheme. In other cases, Lax-Wendroff (LW) or Euler-backward (EB) 
schemes are used to take advantage of the selective dampening of high fre- 
quency modes by these schemes. Although these schemes may be expected to 
reduce the noise generated at the mesh interface, they require a doubling 
of computer time. It is therefore essential that the mesh interface be 
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carefully handled to minimize, or control, the noise near the interface 
and avoid using these expensive time-differencing schemes. If the noise 
does not seriously contaminate the larger scale solution, filtering before 
output will adequately remove the short-wave components. 

The above discussion has primarily been concerned with the handling of 
the inert io-gravity and computational modes inherent in limited-area, pri- 
mitive equation models. The use of primitive-equation models is most 
essential for phenomena in which the Rossby number may approach unity. 
Several groups active in fine-mesh modeling have been testing an implicit 
treatment of the gravity waves present in these models which will permit 
an increase in the time step. The semi-implicit time differencing, after 
Robert (1969), implicitly treats both the spatial gradient of geopotential 
height in the equations of motion, and the divergence term in the continuity 
equation that is multiplied by the mean geopotential height. As this re- 
quires a simultaneous solution for the geopotential and horizontal 
velocity components, the equations may be combined to form a Helmholtz equa- 
tion (usually in (f)) , which may be solved by standard techniques. Solution 
of the Helmholtz equation requires additional computer time, but the 
implicit treatment of the gravity modes permits a much longer time step. 
Several semi-implicit treatments of barotropic models (e.g. Chen and Miya- 
koda, 1974) have documented a savings of a factor of about five over 
explicit treatments. In the baroclinic case a set of coupled Helmholtz 
equations results, which requires special treatment (Robert et al, 1972; 

Sela and Scolnik, 1972). In particular, the convective heat release schemes 
may be troublesome in semi-implicit models. Because of these difficulties, 
and the reduction in phase speed of gravity modes that are the primary 
characteristic of the primitive equation set, some groups have reexamined 
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the use of filtered equation models (e.g. Environmental Prediction Research 
Facility; Bengtsson, 1974). 

The existence of direct-solvers of Poisson or Helmholtz equations 
(work in progress at EPRF and elsewhere) may also contribute to the time 
advantage of the semi- implicit scheme over explicit time -differencing. If 
it could be shown that no loss of accuracy occurs with the semi- implicit 
solution, the time advantage could be exploited to run fine-mesh models 
with existing computer systems. This is one of the most important reasons 
that tropical cyclone models may be feasible at the present time. 

It might be noted that a savings of about one half can be achieved in 
explicit models by arrangement of the solution order to permit averaging 
the pressure gradient terms in the equations of motion (J. Brown and K. 
Campana, NMC technical note). That is, the thermodynamic and j) 1000 (in 
pressure-coordinate models) equations are first solved to calculate the 
new geopotential field. Proper weighting of the pressure gradient terms 
at the past, present and future time steps will evidently permit a doubling 
of the time step. Although some additional storage is required, no Helm- 
holtz equation must be solved, so the full time advantage is realized. 

Arakawa et al (1974) have studied the effect of space-differencing on 
two computational problems that arise in simulating large-scale atmospheric 
motions with the primitive equations. The first problem is the proper 
simulation of geostrophic ad justment-the process by which the energy of 
locally excited inert io-gravity waves is dispersed away into the wider 
space. A second computational problem is the simulation of the slowly 
changing quasi-geostrophic motion after it has been established by geos- 
trophic adjustment. The letters A-E for the space difference schemes 
listed in Table 4 refer to the five orientations of the dependent variables 
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illustrated in Fig. 1. Arakawa et al (1974) suggest that schemes B and C 
have more satisfactory dispersive properties for inertio-gravity waves, and 
should simulate geostrophic adjustment for meteorological waves better than 
the other schemes. A careful treatment of the inertia terms in the 
momentum equation with the grid schemes B and C was shown by Arakawa et al 
(1974) to prevent a computational cascade for purely non-divergent flow. 
Thus this arrangement can be expected to provide a proper representation 
of the large scale, quas i-nondivergent atmospheric motions. 

The majority of the models listed in Table 4 use scheme A, which 
carries geopotential and velocity values at every grid point. Scheme D, 
which has the u(v) velocity points between meridional (zonal) geopotential 
values, has similar properties as Scheme C for preventing cascade, but is 
less satisfactory in terms of geostrophic adjustment. The two research 
models by Chen and Miyakoda (1974) and Bleck and Shapiro (1974) do make 
use of the more favorable Scheme C. Chen and Miyakoda point out that 
Scheme C is also very convenient in terms of applying semi- implicit time 
differencing described above. Elvius and Sundstrom (1973) have also in- 
vestigated the accuracy and stability of the shallow water equations using 
semi- implicit time differencing with Scheme C. In its most compact form, 
the u and v gridpoints alternate locations between odd and even time steps. 
In this case a single Helmholtz equation must be solved, rather than two 
disconnected Helmholtz equations if the velocity components do not alter- 
nate grid points with time steps (Elvius and Sundstrom, 1973). It thus 
appears that Scheme C with semi- implicit time differencing has a number of 
advantages for fine-mesh modeling. 

Various coordinate systems may be used to represent the vertical struc- 
ture of the atmosphere. Kasahara (1974) reviews the features of the 
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primitive equations in various vertical coordinate systems. The pressure 
coordinate system, which is used in most of the models listed in Table 1, 
has certain computational limitations in the vicinity of mountains. Use 
of the vertical height (Z) coordinate was not common until the system was 
adopted by the NCAR general circulation group. The limited-area models by 
Williamson and Browning (1974) and by Kreitzberg et al (1974) use input 
from the NCAR GCM, and thus are cast in height coordinates. Phillips (1957) 
proposed a modified version of the pressure coordinate, the sigma system, 
to make the earth 1 s surface a coordinate surface. This system is used by 
many operational numerical weather prediction centers, and is the coordinate 
system for the LFM. The use of potential temperature (0) as a vertical 
coordinate has advantages for resolution of strong vertical gradients of 
wind. Shapiro and Hastings (1973) have shown that details of frontal 
structure can be resolved in 0 coordinates using the synoptic-scale network 
over the United States. Furthermore, development of the high wind speed 
regions associated with frontogenesis seems to be better predicted in 
0 coordinate models than in pressure coordinate models (Bleck and Shapiro, 
1974). However, the difficulty in handling the lower boundary condition 
in isentropic coordinates, and the lack of adequate resolution of rawinsonde 
stations over the oceans, will probably limit the usefulness of 0 coordinate 
models to continental regions with good data coverage. 

Initialization . The specification of initial conditions is one of the 
most important problems for integration of the primitive-equation models, 
because of the gravity waves that are excited if the mass and wind fields 
are not properly balanced. Although the geostrophic adjustment process will 
eventually disperse these waves, the associated divergence and pressure 
fluctuations may seriously contaminate the meteorological modes during the 
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the transition phase. The precise effect on the large-scale motion fore- 
cast, or on such derived quantities as precipitation amounts, has not been 
evaluated. These effects are probably model-dependent, especially with 
regard to processes that depend on the divergent portion of the motion. 
Recent developments relative to initialization of primitive-equation models 
are discussed by Haltiner and Williams (1974), including the developments 
in four-dimensional assimilation of data. Insertion of new data during the 
forecast also generates inertio-gravity waves in the model. Thus careful 
treatment of data prior to, or during, the integration is required to 
minimize adverse effects on the forecast. 

Many of the fine-mesh models listed in Table 4 follow the practice of 
filtered models, in that the wind fields are derived from a solution of 
the balance equation with an objec ively analyzed geopotential field. As 
shown by Phillips (1960), inertio-gravity waves will be generated in this 
case, unless an appropriate divergent wind component is included. The 
required form of the divergent wind component is far from clear, as several 
operational centers have found that the use of a co-equat ion, which is 
seemingly appropriate to the balance equation, does not improve the initial 
adjustment phase in the model, or the forecast results. For fine-mesh 
models, and particularly for tropical applications, it seems rather clear 
that the traditional initialization methods are inappropriate. It can be 
shown that the mass field tends to adjust to the wind field as the scale of 
the motion is reduced. For mesoscale models the mass field immediately 
tends to adjust to the specified wind field, rather than vice versa (Anthes , 
1974). Furthermore, as the distance between observations is decreased, 
the errors in determining the geopotential approach the magnitude of the 
gradient. Consequently, in mesoscale studies, as in the tropical cases, the 
analysis and initialization phase must consider both wind and mass fields. 
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One must question the mesoscale initialization technique described by 
Rao and Fishman (1974), because the wind field was forced to adjust to 
the mass field. This study does suggest that a staggered grid (Scheme D 
in Fig. 1) improved the convergence rate over the original Miyakoda and 
Moyer (1968) iterative technique. Likewise, Paine and Kaplan (1974) use 
operationally-analyzed height fields interpolated to a 42-km grid, and 
calculate goestrophic winds. The imbalance created by this technique 
leads to massive readjustment through large-amplitude, rapidly propagat- 
ing gravity waves. It might be noted that the tropical model of 
Krishnamurti (1973) is initialized from the wind field, from which the 
geopotential field is derived. The divergent component is derived from 
an omega equation. Williamson and Browning (1973) initialized 2 1/2°, 

1 1/4° and 5/8° versions of the NCAR LAM with data interpolated from a 
5° global model. Even though the original mass and wind fields were 
dynamically consistent, the fields are not necessarily in balance after 
the interpolation process. Even more important is that no meteorologi- 
cally significant features with scales less than 5° are included in the 
data. Normally a 5° grid, or at most a 2.5° grid, would be sufficient 
to forecast these larger scale features. One justification of forecast- 
ing with 1 1/4° and 5/8° grids would be the anticipated development of 
small-scale features. Even in this case, it would appear essential that 
the maximum amount of data is used in each grid, rather than simply in- 
terpolating from larger scale data. 

Temperaton (1973) tested various methods that have been proposed for 
initializing primitive-equation models. Temperaton concludes that forc- 
ing the wind field to adjust to the mass field is "unnatural" for the 
small-scale components. Rather provision needs to be made for the mass 
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field to adjust to the wind field, as well as vice versa. The best method 
of accomplishing the dynamic initialization appeared to be a cycle of 
averaging two forecasts integrated forward and backward in time from an 
initial state, alternately restoring ({) or the wind components during each 
iteration. Temperaton demonstrated that it was not sufficient that the 
initialization method just reduce the root-mean square error for the geo- 
potential field. Rather the principal test of the method must be in 
elimination of the unwanted gravity-wave modes during the integration phase. 
A suitable measure of these waves is the divergent component or vertical 
motion. 

The difficulty in applying these dynamic initialization techniques is 
that the computer time required for the forward and backward integration 
of the model is more or less equivalent to that of the actual forecast. 
Consequently the variational approach originally developed by Sasaki (1958) 
and applied by Lewis (1972) for operational use may have considerable 
promise. With this approach the objectively analyzed wind and temperature 
fields are simultaneously modified, based on dynamical constraints. These 
constraints must be cast in the same form as in the prediction model, but 
the model itself is not integrated. The variational approach requires 
solution of a Helmholtz equation, but with the development of direct 
solvers for this type of equation, as mentioned above, the initialization 
phase need not consume a large fraction of the actual forecast time. 

Although arbitrary weighting factors must be specified to determine the 
relative importance of observations to dynamical constraints, these factors 
may be used to emphasize the wind field in the desired regions. 

Representation of physical processes . One of the prime motivations 
for fine-mesh modelling is a more accurate representation of the physical 
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processes. These may involve flow over detailed orography, momentum and 
energy redistribution in convective systems, or the vertical flux of 
properties through the planetary boundary layer to the upper atmosphere* 
Methods of representing these processes are presently being actively re- 
searched. This is particularly true for the effects of convective-scale 
heat release on the larger scale flow and for boundary layer processes. As 
these processes occur on scales smaller than the grid size, their effect 
on the grid scale flow must be parameterized. The success of these 
parameterization schemes will be increased if the large scale processes 
are represented on the scale that triggers the exchange, and/or receives 
the primary feedback signal. This will require fine-mesh modeling as the 
triggering processes appear on scales no larger than the mesoscale. 

Garstang and Betts (1974) has summarized the extensive research regard- 
ing representation of latent heat release in both prognostic and diagnostic 
studies. One of the key features is an improved understanding of how the 

large scale responds to the presence of a cloud, in particular through 

the subsidence heating in the environment in response to the buoyant ascent 
of cloud air. However, the vertical distribution of this heating is 
dependent on the entrainment /detr a inmen t process, and the associated down- 
draft-production mechanism, which are not well understood. The very active 
research in this area, and the results of further studies based on GATE 
data obtained during 1974, promise an improved parameterization of the latent 
heat release. The basic framework for these schemes has been established 
by studies such as Ooyama (1971) and Arakawa and Schubert (1974). In these 

studies the two key model elements are the large scale forcing of the 

cumulus clouds and the feedback to the environment through a spectrum of 
cloud sizes. Arakawa and Schubert (1974) particularly emphasize the role 
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of the boundary layer processes in determining the amount of mass ascent 
through clouds. An example of a fully coupled boundary layer and latent 
heat parameterization model has been described by Arakawa and Randall 
(1974). Kreitzberg et al (1974) also anticipate a dominant role for the 
planetary boundary layer model, with the coupling of the horizontal flow 
to the detailed vertical structure requiring a mesoscale rather than 
cyclonic scale resolution. Whereas Arakawa and Randall use a simple 
boundary layer model primarily to improve the vertical flux calculation 
at the lowest levels in the global model, Kreitzberg et al actually insert 
about 15 additional levels and calculate the vertical flux at each level. 

In general the models described in Table 4 incorporate rather crude 
schemes for parameterization of latent heat release. Those derived from 
operational models use a form of convective adjustment to remove insta- 
bility. Such schemes cannot be directly used in regions that are typically 
conditionally unstable--such as the tropics, in the warm sectors of extra- 
tropical cyclones, or in continental regions during the summer. If the 
convective adjustment scheme is not modified, the conditionally unstable 
regions will be completely eliminated, with a large overestimate of preci- 
pitation. One model that includes an extremely complex latent heat para- 
meterization scheme is Kreitzberg et al (1974). The potential instability 
at each gridpoint is released by a detailed, one-dimensional cloud model 
that includes microphysical processes. This brute-force technique is 
extremely costly in terms of computer resources, but may be particularly 
useful if the simulation or modification of cloud processes is being 
studied. An intermediate example is the use of a modified Kuo parameteri- 
zation scheme by Krishnamurti et al (1973). This scheme has been tested 
in various tropical systems and in some mid-latitude diagnostic studies. 
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It appears that the association of large-scale and boundary layer conver- 
gence of moisture with convective mass flux gives useful results. Even 
though the physical interpretation of the vertical heat distribution 
mechanism is questionable for the Kuo scheme, the utility of the scheme 
for short-term forecasts of convectively-active weather systems has been 
established. 
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Fig. 1 Various spatial arrangements of variables suggested by Arakawa 
et al (1974), as described in the text and Tables 1 and 4. 
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Fig. 2 An example of a meshed grid with a discrete change in resolution 
by a factor of two. 
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Fig. 3 An example of a triply nested grid arrangement (Harrison, 1973). 
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